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I. BACKGROUND 

The Keldysh nonequffibrium Green function technique is used very widely to describe transport phenomena in 
mesoscopic systems. The technique is somewhat subtle, and a rigorous treatment would require much more timeTha^ 
we have at our disposal, see, for example, the text-book by Haug and Jauho [1]. The goal of these lectures L to 
gve a rough idea of how the technique works, and illustrate it with a few simple examples. Our main S 
iS*} 6 T cun f « diode; this calculation can be done in a number of dEnt 

ways, such as the Landauer-Butt.ker formalism, or by a Master Equation approach, and our hope is that studymg 
this exactly solvable problem by different methods will give insight to the different formalisms. 

II. REVIEW OF PERTURBATION EXPANSION OF THE GREEN FUNCTION 

rLtt pZ£! 5^ tbe a ^ n M iS f r mar ^, equilibrium Green 011 the kwl *s discussed in Many- 

wl l^TS I (Transport m ^structures) [2]: Nevertheless, I summarize some of the central features here. 
We start by defining a time-ordered (also called causal) zero-temperature singlelparticle Green function 

n (* |*o> ' W 

Next, recall the various representation pictures of quantum mechanics: 

. Schrodinger picture: the wavefunctions are time-dependent: i^{t) = H^t); the operators are constant. 
' *%?™Q ) *% t pictmei ^ wavefunctions are constant; the operators are time-dependent: 0(t) = 

" f tBI * C tion P ictl ^ e: ^functions develop under the influence of the "difficult" interaction part of the 
Hamiltonian if = H + V: #t) = e**°V^(0). The operators develop under the influence of the W 
non-mteracting Hamiltonian H only: <D(t) = e iHot 6(Q)e~ iHat . [9] 

The time-development of the interaction picture is often expressed in terms of a U operator (or matrix): 

i>(t) = U(t)i,(0), «7(t) = e iff ° t e -^ t . (2) 
We next introduce the 5-matrix, which changes the wavefunction from i>{t') to tp(t): 

#) = 5(i,^(0- (3) 

Thus, 

S(t,t) = U(t)tf(t'). (4) 

The 5-matrix obeys the group property S(t,t') = S(t,t")S(f,f). Itothermore, the 5-matrix can be expressed as a 
time-ordered product [for a derivation based on the equation of motion for 5, see [2],[3],[4]: jessed as a 

5(t,i , ) = Te- i /«' d *i^). (5) 

the Gree * * nM * iQa contai f a diffiai%: » involves the exact ground state of the system. But this 
is one of the hmgs we want to compute with the Green function! To make any progress, we must express the™ 

Seltt^ 

|*o) = 5(0,-oo)|$ >, ( 6) 
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The proof for this relation is rather subtle (see Ref.[5], pp. 61-64) and we do not reproduce it here. In analogy with 
(6), we have 

<*o| = <$o|5(co,0). (7) 

The object ($015(00, 0)5(0, -oo)|$ ) may have an ill-defined phase, this is, however, canceled by a similar phase-factor 
arising from the numerator of the definition of the time-ordered Green function: 

G{*,t; X ;t>) = ^^oir{g(^-°o)^(«,«)f t y,*')}|go) . (8) 

($015(00, -co) |$ > ; 

This important result generates the systematic perturbation scheme for the Green function. The calculation proceeds 
by expanding the 5-matrix (both in the numerator and in the denominator) in V(t): 

5(oo, -oo) = £ i-ip- / dty dt n T{V(h) ■ ■ ■ V(tn)} , (9) 

n=0 Ul J -°° 

Since each V contains three or four field operators, we need a device for evaluating expectation values such as 

(SolT^W^^^O^fe^fe)^^)}!^} , (10) 

and higher order similar terms. These expressions are evaluated with Wick's theorem, which says that the result of 
(10) is the sum of all pairwise contractions. Thus (10) gives rise to six terms: for example, ij>(t) can be paired with 
$(t'), and the remaining 4 operators can be paired in two different ways. Thus one gets 3x2 = 6 terms. These 
terms are most easily expressed in terms of Eeynman diagrams. 

The various diagrams have quite distinct properties. Some of them are disconnected diagrams. These are ex- 
actly canceled by the denominator in (8). This is good because disconnected diagrams often diverge! After some 
combinatorics (see Sect. 12.3.3 in [2]), one finds that 



-+0O 

G{x,f,x\t') = -*£(-*)" / *l • ■■dt n {$ \Ti,(x,t)^(x',t)V(t 1 ) • • -F(« n )|$ >conn , (11) 



where the summation only includes topologically different connected diagrams. (11) is the desired perturbation 
expansion for the Green function, and it forms the starting point for many calculations. 



III. CONTOUR ORDERED GREEN FUNCTIONS 

The central quantity in constructing the perturbation theory for Green functions is the 5-matrix 5(oo,-oo). In 
nonequilibrium there is no guarantee that the system returns to its initial state for asymptotically large times. In 
fact, often it does not. Consider, for example, an important problem in surface physics, where atoms or molecules 
impinging on a surface exchange charge with the surface, and hence the initial state at i = — oo is very different 
from the final state at t = +00. Thus, one should avoid any reference to the asymptotically large times in the non- 
equilibrium theory. The general formulation of the theory is slightly more complicated than in the equilibrium case. 
As we shall see, however, the abstract structure of the theory bears a close resemblance to the equilibrium theory. 

The nonequilibrium problem is formulated as follows. We consider a system evolving under the Hamiltonian 

H = h + H'(t). (12) 

Here the time-independent part of the Hamiltonian h is split in two parts: h = H + Hi, where H is "simple" (in the 
sense that it can be diagonalized,- and hence Wick's theorem applies) and Hi is "complicated" (in the sense that it 
contains the many-body aspects of the problem, and hence requires a special treatment). It is further assumed that 
the nonequilibrium part vanishes for times t < to. The nonequilibrium part could be, e.g., an electric field, an light 
excitation pulse, or a coupling to contacts at differing (electro) chemical potentials. 

One often lets to -* -00 at a suitable point. This procedure simplifies the treatment, and in order to display the 
structure of the nonequilibrium theory as concisely as possible, we take this limit. However, in doing so the discussion 
of transient phenomena is excluded, and must be treated separately. 
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Before the perturbation is turned on, the system is described by the thermal equilibrium density matrix, 

6{h) = Ir[exp(-/4 • (13) 

The task is to calculate the expectation value of a given observable, to which one associates a quantum mechanical 
operator O, for times t > to: 

(0(t)) = Tr[ Q (h)0 H {t)}. (14) 

The subscript H indicates that the time-dependence is governed by the full Hamiltonian, i.e., O is written in the 
Heisenberg picture. The definition (14) can be generalized to two-time (or n-time) quantities (Green functions 
correlation functions) in an obvious fashion. . ' 

One should note that we use an equilibrium density matrix g(h) in (14), and not some time-dependent q. Physically 
this means that the thermodynamic degrees of freedom, contained in h, do not follow instantaneously the rapid 
variations contained in H'(t). Other choices may be possible, but we make this choice here because of the difficulties 
related to other choices [see the discussion on pp. 214-216 in[4]. 

The general plan of attack is similar to the equilibrium case. We transform the "hopelessly complicated" time- 
dependence of O h to a simpler form, namely to that of So . Since there are two operators to be eUminated, i.e., 
the time-dependent external perturbation H'(t), and the "complicated" interaction term H u we expect to meet more 
complicated transformations than in the equilibrium case. However, with suitable generalizations, it can be shown 
that the nonequilibrium and equilibrium formalisms can be made structurally equivalent. 

We skip here the technical details of these transformations because the final result is intuitive; an interested 
reader*may consult [1]. Instead of the time-ordering operator used in equihbrium theory, we now introduce the 
contour-ordering operator, which orders the time-labels according to their order on the (Keldysh) contour. With this 
generalization, the topological structure of the perturbation theory is identical to the equilibrium theory. Thus the 
contour-ordered Green function, defined as ' 

0(1,10 & -i(Tc v bPitm H (l')}) , (15) 
satisfies the same Dyson equation as the equilibrium function: 

G(l, 10 = <?o(l, 10 + J d 3 x 2 jf dr 2 G (l, 2)t/(2)G(2, 10 + J dh 2 J d3x 3 J dr 2 J dr 3 G (l, 2)S(2, 3)G(3, 10 , 

where all the time-labels now reside on the contour shown in Fig!. In writing (16), we assume that the nonequilibrium 
term in the Hamiltonian can be represented by a one-body external potential U. The interactions are contained in 
the (irreducible) self-energy S[G]. The contour representation is rather impractical in calculations, and we want to 
replace the contour by real time integrals. This procedure is called the analytic continuation. 

As mentioned above, a simplification occurs if we can set t -co. If the interactions are coupled adiabatically, 
the contribution from the [t ,t - i0\ piece vanishes. The information lost by this procedure is related to initial 
correlations. In many physical situations, for example, in steady state transport, it appears plausible that the initial 
correlations have been washed out by the interactions when one reaches the steady state. On the contrary, if one 
studies transient response, the role of initial correlations can be important. This represents an interesting, and rather 
open problem. Here we consider the t -> -co limit. The contour thus consists of two parts, from -co to co, and 
from oo to -co. Any time residing on the first part (call this Q) is earlier in the contour sense to any time residing 
on the latter part C 2 . The contour-ordered Green function thus contains four different functions: 

fG e (l,10 ti,t v eCi 

Gfl l'l = I G> ^' *i e °»> 'i' G C i „ « 

{G e {i,V) t u t v ec 2 
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Here we have introduced the causal, or time-ordered Green function G c , 

G c (l,l') = -i{T[fe(l)^(l')]) (18) 
= -ie{h - t y )ty H (l)tf H (l>)) +i9(t v -ii)(^(l')V/f(l)) , 



the "greater" function G > , 
the "lesser" function G < , 



G>(l,l') = -i{V'fi(l)^(l')), (19) 



G<(l,l') = +i(^(l')^(l)), (20) 
and the antitime-ordered Green function Gg, 

G c -(1,1') = -i(T[fe(l)^(l')]> ' (21) 

= -ie{t v - h){i> H {i)^ B {\')) + iefa - i!<)ft4(i>iMi)> • 

Since G C + G 5 = G < + G > , there are only three linearly independent functions. This freedom of choice reflects itself 
m the literature, where a number of different conventions can be found. For our purposes the most suitable functions 
are the functions G><< (which are often also denoted by a common name, "correlation function"), and the advanced 
and retarded functions defined as 

G°(l,l') = tf(fi,-< 1 ){WH(l) I ^(l')}> 

= 8{tv - h)[G<(l, 1') - G>(1, 1')] , (22) 

and 

C(i, i') = -aft - M<ftMi),^(i')}> 

= 0(h - t v )[G>(l, 1') - G<(1, 1')] . (23) 

Here early brackets denote an anticommutator. We observe that G r - G a = G > - G K . 

The main technical trick needed in the analytic continuation is included in the following theorem. 
Langreth Theorem. In considering the Dyson equation (16) we encounter terms with the structure C = AB or 
explicitly, ' ' 



C{t u tv)= j drA(t lt T)B(T,tv) 
Jc 



(24) 



and their generalizations involving products of three (or more) terms. Since we are presently only concerned with 
temporal variables, we suppress all other variables (spatial, spin, etc.), which have an obvious matrix structure. To 
evaluate (24) let us assume for defmiteness that ti is on the first half of the contour ("on the way out") and that t v is 
on the latter half ("on the way back"). In view of our discussion in connection with (18)-(23), we are'thus analyzine 
a "lesser" function. J B 

The next step is to deform the contour as indicated in Fig. 2. Thus (24) becomes 

C^fr.tiO = J c dTA(t lt r)B<(T,tv) + dT J 4<(t 1 ,r)B(r,t 1 ,) . ( 2 5) 

Here, in appending the sign < to the function B in the first term, we made use of the fact that as long as the 
integration variable r is confined on the contour Gj it is less than (in the contour sense) t v . A similar argument 
applies to the second term. Now consider the first term in (25), and split the integration 'into two parts: 

f drAit^B^ty) = P dtA>(t 1 ,t)B<(t,t 1 ,)+ r°°dtA<(t 1 ,t)B<(t,t 1 A 

/CO 
dtA r {t u t)B<{t,t v ), (26) 
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FIG. 2: Deformation of contour C 



where we used the definition of the retarded function (23). A similar analysis can be applied to the second term 
involving contour C{; this time the advanced function is generated. Putting the two terms together, we have the first 
of Langreth's results: 



/oo 
dl{A r {t l ,^B<(t,t s .) + A < (t lt t)B a ( L t,t 1 .)] 
CO 



(27) 



The same result applies for the "greater" function: one just replaces all <'s by >'s. 

It is easy to generalize the result (27) for a (matrix) product of three functions: If D = ABC on the contour, then, 
on the real axis, one has 

D< = A T B r C< + A T B<C a + A<B a C a . (28) 

Once again a similar equation holds for the "greater" functions. 

One often needs the retarded (or advanced) component of a product of functions defined on the contour. The 
required expression is derived by repeated use of the definitions (18)-(23), and the result (27): 

/CO 
dt[A r (B> - B<) + (A> - A<)B a ] 
■OO 

= 6(h - t v ) [y* 1 dt(A> - A<)(B> - B<) 
+ f 1 ' dt(A> - A<)(B< - B>) 

J —CO 

= r dtA^tu^B^tv) . (29) 

In our compact notation this relation is expressed as C = A r B r . 

When considering the various terms in the diagrammatic perturbation series one may also encounter terms where 
two Green function lines run (anti)parallel [this is the case, for example, of a polarization or electron-phonon self- 
energy (parallel Fermion and Boson Green functions)]. In this case one needs the "lesser" and/or retarded/advanced 
components of structures like 

C(t,t') = ^rjr'^T.r'), 

D(t,t') = A(r,r')B(T',T), ( 30 ) 

where r and t' are contour variables. The derivation of the required formulae is similar to the analysis presented 
above, and is left for an exercise. One finds 

C<(t,f) = A<(t,t')B<(t,t') , 

D<(t, = A<(t,t')B>(l/,t) , ( 31 ) 
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TABLE I: Rules for analytic continuation 



Contour Real axis 



C = f c AB 


C < = S t [A r B<+A < B a ] 




C r =f t A r B T 


D = f G ABC 


D< = J t [A r B r C < + A r B < C a + A < B a C a ] 




D r = f t A r B r C 


C(t,t') = A(t,t')B(t,t') 


C<{t,H) = 4<(t,0B < (t,f) 




C-(t, f ) = A<(t, t')B r (t, t') + A r (t, i/)B<(t, f) 




+A r (t,t')B r (t,t') 


D(t,t')=A(t,t')B(t',t) 


D<(t,i') = A<(t,-e)B>(i/,t) 




D r (t,t') = A<(t^)B a ^,t) + A r (t,t')B<{f,t) 



and 

C r (t,t') = A < (M05W) + ^W)B < (M0 + 4 r fti0£''(i,i')> 

D r (t,t') = A r (t,t')B<(t,t)+A<(t,t')B a (t',t) 

= A<(t,t')B a (t',t)+A r (t,t')B<(t',t) . (32) 

As earlier, the relations (31) can immediately be generalized to "greater" functions. For a quick reference, we have 
collected the rules provided by the Langreth theorem in Table 4.1. 



Equilibrium Electron-Phonon Self-Energy. 

The retarded electron-phonon self-energy £ ph is a central object in the analysis of many physical properties of metals and 
semiconductors. At finite temperatures one conventionally uses the Matsubara technique to perform the analytic continuation 
[for an extended discussion, see [4]. The Langreth theorem can be used to give a very compact derivation of Sp h . In lowest 
order in the electron-phonon matrix element M q we have 

S P h(ft,r,r') =iY,\M<,\ 2 G(k-q,TS)D{ q ,T,T') . (33) 
t 

Here G is the free-electron Green function while D is the free-phonon Green function. Equation (33) is in a form where we 
can apply (32). In equilibrium all quantities depend on time only through the difference of the two time-labels, and it is 
advantageous to work in frequency space. Performing the Fourier transform gives 

W") = i / :§ Y, KI 2 [ G< (*-ff.*> -e)D r -(q,e)+ G r (k - g,u - e)D < (g,e) + G r (k-g,u - e)D r (q,e)] . (34) 

The expressions for the free equilibrium Green functions are (the reader is urged to verify these relations!): 

D < (g,u) = -2m[{Nq + l)S(u+u q ) + N q 5{u-u q )], 
D r (g,w) = 1 : 1 , 

W - U>q + IT) li)+lil q + irj 

G < (k,w) = 2mnF(u>)6(u> — Eh) , 
Substituting these expressions in (34), one finds after some straightforward algebra 



Sph(fc^) = YX -"'(«*-») + 1 + + "*(■*»-«) 

[W — U> 4 — Sk-q + iri W + Uq- Sk-q + it] . 



(36) 



The shortness of this derivation, as compared to the standard one, nicely illustrates the formal power embedded in the Langreth 
theorem. 
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IV. KBLDYSH FORMULATION 

We now derive the Keldysh integral equation for the correlation function G<, by applying the analytic continuation 
rules on the complex time Dyson equation. Thus, applying the rule (28) on (16) [10] yields 

G < = G$ + G r 2 r G < + G r X<G a + GZVG a . (37) 
We proceed by iteration with respect to G<. Iterating once, and regrouping the terms we obtain 

G< = (1 + GSSHGtfa + S°G a ) + (GJ + G r X r G r )E<G a 

+ G$rGV?G< . (3 8 ) 

The form of (38) is very suggestive, and it easy to convince oneself that the infinite order iterate is 

G< = (1 + Grs r )G$ (1 + £°G a ) + G r S<G a . (39) 

One can show, using the Dyson equation for the retarded Green function, that the first term in (39) vanishes for 
steady state systems, if the system was in a noninteracting state in infinite past. Thus, in many applications (also as 
is done below) it is sufficient to only keep the second term. Despite of its simple appearance, one should note that the 
Keldysh equation is an integral equation for the correlation function G<, because usually £< is a functional of G<. 



V. NONEQUILIBRIUM TECHNIQUES IN MESOSCOPIC TUNNELING STRUCTURES 

Let us now examine how the techniques discussed above can be used to derive a perturbation scheme for mesoscopic 
systems far from equilibrium. We once again recall that the basic difference between construction of equilibrium and 
nonequilibnum perturbation schemes is that in nonequilibrium one cannot assume that the system returns to its 
ground state (or a thermodynamic equifibrium state at finite temperatures) as t -» +00. Irreversible effects break 
the symmetry between * = -00 and t = +oo, and this symmetry was heavily exploited in the derivation of the 
eqmhbnum perturbation expansion. In nonequihbrium situations one can circumvent this problem by allowing the 
system to evolve from -00 to the moment of interest (for definiteness, let us call this instant to), and then continues 
the time evolvement from t = t Q back to t = -00. The advantage of this procedure is that all expectation values are 
define* with respect to a well-defined state, i.e., the state in which the system was prepared in the remote past. The 
price is that one must treat the two time branches on an equal footing. 

In the context of tunneling problems the nonequilibrium formalism works as follows. In the remote past the contacts 
(i.e., the left and right lead) and the central region are decoupled, and each region is in thermal equilibrium. The equi- 
librium distribution functions for the three regions are characterized by their respective chemical potentials; these do 
not have to coincide nor are the differences between the chemical potentials necessarily small. The couplings between 
the different regions are then established and treated as perturbations via the standard techniques of perturbation 
theory, albeit on the two-branch time contour. It is important to notice that the couplings do not have to be small 
e.g., with respect level to spacings or k B T, and typically must be treated to all orders. 



VI. MODEL HAMILTONIAN 

We split the total Hamiltonian in three pieces: H — H c + Ht + H^, where H c describes the contacts, H T is 
the tunnehng coupling between contacts and the interacting region, and H cea models the interacting central region, 
respectively. Below we discuss each of these terms. 

Guided by the typical experimental geometry in which the leads rapidly broaden into metallic contacts we view 
electrons m the leads as non-interacting except for an overall self-consistent potential. Thus, the contact Hamiltonian 
is 

H o= Yl e *<* c L c *a, (40) 
k,aeL,R 

and the Green functions in the leads for the uncoupled system are: 

ti a (t ~ *0 = i(4a(t')^ a (t)) = if(e° ka ) exp [ - is ka (t - 1')] , ( 41 ) 
flta(*~«0 = ^±*=F*0{{c to (*),cl o (t')}> = T itf(±t T «0«ap[-ie»a(«-*')]- (42) 
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Here f(e ka ) - [exp[(e ka - n a )/k B T\ + 1] is the equilibrium distribution in a given lead. 

The coupling constants between the leads and the central (interacting) region depend, in principle, on the actual 
charge densities (accumulation and depletion regions, charge build-up in the central region, etc.), and they should be 
determined via a self-consistent calculation. We do not attempt to perform such an analysis, and assume that these 
parameters are known. Thus we write 

ffr= X) IWL4 + 1I.C.]. (43) 

n 

Here, {d£ } and {d n } are the single-electron creation and aimihilation operators for the complete and orthonormal 
set of the states \n) in the interacting region. 

The form chosen for the central region Hamiltonian depends on geometry and on the physical behavior being 
mvestigated. We will discuss three particular examples in detail. In the first, the central region is taken to consist of 
non-interacting levels, 

•Hcen = ^Emd^dm . (44) 
m 

Here, <4(d m ) creates (destroys) an electron in state m. The choice (44) represents a simple model for resonant 
tunneling. Below we shall present general results for an arbitrary number of levels, and analyze the case of a single 
level, which is exactly solvable in detail. 
A second example of interest is resonant tunneling with electron-phonon interaction: 

= £od f d + dtd M g i & l+ a_ f ] . (45) 

Here, the first term represents a single site, while the second term represents interaction of an electron on the site 
with phonons: a^af) creates (destroys) a phonon in mode q, and M ? is the interaction matrix element. The Ml 
Hamiltonian of the system must also include the free-phonon contribution H pb = V-Tja^-ataj. This example while 
not exactly solvable, is helpful to show how interactions influence the current. This model is presently under intense 
study because of its relevance to inelastic effects in molecular electronics. A simple discussion is available in Refill. 
Thecal example consists of an Anderson-type model for electron-electron interactions in the central region: 

■Hcen = ]T) £odJd ff + Un^ni . (46) 

Here a is a spin-label, n„ is the occupation number operator of spin-state tr, and U describes the on-site Coulombic 
repulsion. This model has been the topic of intense study because of the very rich (and complicated) physics that it 
describes. For example, in the low temperature limit it exhibits Kondo behavior. In the high temperature limit the 
equation-of-motion technique allows a relatively simple analysis, see [1]. 

VII. GENERAL EXPRESSION FOR THE CURRENT 

The current from the left contact through left barrier to the central region can be calculated from the time evolution 
of the occupation number operator of the left contact: 

JL = -e{N L ) = -^{[H,N L }), (47) 

where N L = £ fc ,ae£ c\ a c ka and H = H a + H T + ffcen. Since H c and H^ D commute with N L , one readily finds: 

Jl = J 12 PWfcLA) - V k * a>n {dl,c ka )] . (48) 



Now define two new Green functions (we set ft = 1, and reintroduce it in the final expression for the current): 

G^ o (i-0 = *(ck(tf)d»(t)> , 

s i(4(f)c te (t)>. (49) 



We see that the current is given by the time-diagonal components of the Green functions defined in (49). These 
functions have the property t) = -[G< ka (t,t)]* , and inserting the time-labels, the current can be expressed 

OS 

^ = f M E V ^nG< ka {t,t)\ . (50) 

Next, one needs an expression for G< M (t- 1'). For the present case, with non-interacting leads, a general relation for 
the contour-ordered Green function G n , ka {r,r') can be derived rather easily with the equation-of-motion technique 
Since the nonequihbnum theory is structurally equivalent to equilibrium theory, it is sufficient to consider the T = 
equation-of-motion for the time-ordered Green function G* ka : 

Q 

-^GW* - = **GU(* - *) + £ G »m(* " m£ a , m , (51) 

m 

where we defined the central region time-ordered Green function function Gl m (t - f ) = -ilTfdt (t')d (t)\\ Note 
that it is crucial that the leads be non-interacting; had we allowed interactions in the leads, the ™qua!tion-oimotion 
technique would have generated higher order Green functions in (51), and we would not have a closed set of equations. 

We can interpret the factors multiplying G* nM (t - f ) as the inverse of the contact Green function operator, and 
introduce a short-hand notation: C? nM g^ = By operating with g\ a from right, we arrive at 

GU(* - = £ / dt &m(t - tiW^gUt, - 1') . (52) 
m J 

Since the equilibrium and nonequilibrium theories are topological^ equivalent, this equation has in nonequilibrium 
precisely the same form, except that the intermediate time integration runs on the complex contour: 

G„ m (t, t') = £ / dTiG nn {T, njV^g^in, r') . ( 53 ) 

771 

Here G nm (r, n) is the contour-ordered Green function for the central region. The analytic continuation rules of Table 
1 can now be applied, and we find 

Gn M {t ~ = E/ 4 '^mO - ~ *') + G< m (t - tOflfcto - f')] , (54) 

where the Green functions g<< a for the leads are defined in (41),(42). The Fourier transform of (54) is 

GUM = £ V k \ tm {Gl m {s)gUs) + G< m {e)gUe)\ , (55) 
whereby the current (50) becomes 

Jl = tJ ItM E v ^nV k \ m [G r nm (e)gZ a ( e ) + G< m {e)gt a {e)} } . ( 56 ) 

&nctSn- UBCtUre K iS *° C ° nVert ^ m ° mentum summation s to energy integration and define a level-width 

£fc) . (57) 

where Ba {e) is the density of states. There are two terms in the current expression (56). Consider, for example the 
piece involving GL,„., which we evaluate as fill F ' 



piece involving G r mn , which we evaluate as [11] 
2e t de 

n 

_2e f ds 



/£/ ds " rL (^MG r (^8(e-e k )f L ( e )] = || ^f L (e)r L (e)Re[iG^(e)} 

n J aF 1 *^*)^*) - G °( £ )] • (58) 
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Similar manipulations can be applied to the other term, and the current from left (right) contact to central region 
becomes 

Jlw = 7? / 1 »(r £(fl) (e){G<(e) + f HR) (e) [G'(e) - G°(e)] }) . (59 ) 

Here, the boldface notation indicates that the level-width function T and the central-region Green functions G<< r are 
matrices in the central-region indices m,n. In steady state, the current will be uniform, so that J=J L = -J R and 
one can symmetrize the current: J = (J L + J L )/2 = (J L - J R )/2. Using (59) leads to the general expression for the 
d.c.-current: 

3 = 1 / iH£ r ^ - r *( £ )] G< (s) + [fU^') - Me)r R (e)] [GT(e) - G»(e)]} . (60) 

Often the energy-dependence of the level-width function is not very important, and further simplification can be 
achieved by makuig assumptions on this energy-dependence. In particular, if the left and right line-width functions 
are proportional to each other, i.e., T L (e) = XT R (e), a very sjmple final result can be be achieved. We observe that 
the current can be written as J = xJ L - (1 - x)J R , which gives, using (60): 

3 = ?/ IK rR(e) { [xx ~ (i ~ x)]G<{£) + [xxh - (i - x)f «\ t cr(£ ) - }) • ^ 

The arbitrary parameter x is now fixed so that the first term vanishes, i.e., x = 1/(1 + A), which results in 

^■ ^i^ ^W-^W]}- (62) 

The ratio is weU-defined because the T-matrices are proportional. The difference between the retarded and advanced 
Green functions is essentially the density of states. Despite the apparent similarity of (62) to the Landauer formula 
it is important to bear in mind that, in general, there is no immediate connection between the quantity Tie) and 
the transmission coefficient T(e). In particular, when inelastic scattering is present, there is no such connection. In 
later sections, where we analyze a non-interacting central region, a connection with the transmission coefficient can 
be established. Further, it is possible to derive an analogous result for the average of the time-dependent current [1] 
The results for the current derived in this section [(60) or (62)] form a starting point for many calculations. As 
such, they are powerful formal results, however, one must bear in mind that they involve both the full correlation 
function and the retarded Green functions in the presence of tunneling, and it is by no means obvious that these 
Green functions are easily calculable. The usefulness of (60),(62) depends on whether one is able to devise suitable 
calculation schemes for the Green functions. Below we shall see some examples of how this is done in practice. 



VIII. NON-INTERACTING RESONANT-LEVEL MODEL 



In the non-interacting case, the Hamiltonian is H = H e +H T +H cm , where H cea = £ n e „dt d n . It is straightforward 
to derive the Dyson equation for the retarded Green function: 

G r (e) = g-( £ ) + g^S^G^) , (63) 

where 

*w(e)= 1] Vg a<n g r ka (e)V ka , n , , ( 6 4) 

kaeL,n 

and g r ka is given by (41). Since the lesser self-energy follows explicitly from (64) by analytic continuation, one can 
immediately write down an expression for the lesser function G< using the Keldysh equation (39) for G<: 

G<( £ ) = G'-( £ )S<(e)G»(£). (65) 

In general, the Dyson equation for the retarded Green function requires matrix inversion. In the case of a single level 
the scalar equations can be readily solved. The retarded (advanced) self-energy is ' 

Sr '° (£)= £ FzK^A^flXs), (66) 
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where the real and imaginary parts contain "left" and "right" contributions: A(e) = A L (e) + A R (e) and T(s) 
T L (s) + T R (s). The lesser self-energy is W K ' K ' 



S < (£)= £ \V*?&(e)=i\r*{e)Me)+l*ie)f R (e)]. (67) 



Using the identity 



(G a-1 - G*- 1 ) r(e) ' ^ 68 ) 
where A(e) = i[G r (e) - G a (e)] is the spectral function, one can write G< in a "pseudoequiUbrium" form: 

G<(e) = iA(e)f(e), (69) 



where 



r(e) 

A{£) = le-eo-A(£)]2 + [r( e )/2]2 ' ( 70 ) 
With these expressions the evaluation of the current (60) is straightforward: 

3 = -i/S^ £ ){[ ri ( e )- rK W3/>)-[^(^ L ( £ )-M e )r B ( e )]} 
= t f & r £ (e)r fl (g) 

ft 7 2tt [e - eo - A(e)f» + pT( £ ) /2p [/i(s) ~ /b{e)] * ( 71 ) 

The factor multiplying the difference of the Fermi functions is the elastic transmission coefficient T(e). (Can you 
prove this?) It is important to understand the difference between this result and the result obtained in the previous 
section (despite the similarity of appearance): There, (62) gives the current for a fully interacting system, and 
the evaluation of the retarded and advanced Green functions requires a consideration of interactions (e.g., electron- 
electron, electron-phonon, and spin-flip) in addition to tunneling back and forth to the contacts. To further emphasize 
the difference, let us now suppose that the Green function for the interacting central region can be solved- G r < a (e) = 
[e - e - \{e) ± vy{e)/2) \ where X and 7 /2 are the real and imaginary parts of the self-energy (including interactions 
and tunneling). Then the interacting result for proportionate coupling (62) becomes 

This result coincides with the non-interacting current expression (71) if X(e) -* A(e) and 7(e) -* r(e) = r R (e)+r i (e). 

IX. MASTER EQUATION APPROACH FOR RESONANT TUNNELLING 

The Master Equation (ME) approach is often a quick (and dirty) way of getting to results. In the present, exactly 
solvable case the results of the Keldysh approach and the ME, should of course agree. Let us verify this In ME 
on operates with the occupation probabilities of a given state, P n . In the present case, the two possible states are 
an empty state and "a filled state", P and Pi, respectively. Naturally, P + Pi = 1. To form a link to the Green 
function formulation, we note 

Pi(t) = (*mt)) = ~iG<(t, t) = j ~G<(e). (73) 

Using (69), we calculate as follows (from this on, we work in a so-called Wide-Band Limit (WBL), which means that 
(1) the real part of the self-energy is set to zero, and (ii) the imaginary part, i.e., the T's, are assumed to be energy 
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independent): 



J 2tt (e- eo )2 + (r/ 2 )2 x \ T J 



I 



f • ( 74 ) 

Inthe second line the limit means that we are concerned with the weak coupling limit; this is necessary because the 
ME is inherently low-order perturbation approach (Fermi Golden Rule). Nevertheless, (??) is a bit intriguing result 
because it does not depend on the absolute value of the couplings, only on their relative strength 

Now let us the same using ME approach. Now we simply write down rate-equations for the occupations, using 
"common sense": y ' 5 

Pi = IWPo-rj-.oPi 
P = T^oPi - IWPo 

(75) 

where (note that Golden Rule based T's in the ME are rates, not energies as in the GP formulation) 

MW = T L f L + T R f R 

KT M = T i (l - f L ) + T R (1 - f R ). (76) 

The reader is urged to make sure that she understands the reasoning behind these rates! In the stationary case 
Po = Pi = 0, one readily finds ' 

r = r »-i = r L h+r R f R 
r _i + ri_ r^ + r* 

> r L +v R ' ( 77 > 

in agreement with Eq.(74). 
We can also evaluate the current. In the GF language we had 



de r L r R 

2k (e - e ) 2 + (r/2) 2 f/i(e) " /jj(e)] 
^ p [/i( e o) - /ii(eo)] (78) 



In ME, we calculate 



I = elPoT^-P^} 

= if{F L {l - h) +T R (1 - f R )]I*f L - Q*f L + T R f R )T L (l - h)}} 



e T R T L „ , e rRrL 

= »— f(l - fn)h - f R (l - f L )] = ~~{h - f R ), (79) 

in agreement with the GF result. (Again, do you understand the first line in this calculation') 
EmaUy, it may be of interest to inquire whether one can derive the ME from the GF formalism. The answer should 

be yes , but I m not aware of any elegant demonstrations. The argument could go roughly as given below, but there 

are several places where yellow warning lamps are blinking. You are encouraged to reflect upon them! I leave for you 

the pleasure of finding them. 
For completeness, I give quite a few details. The hamiltonian is 

H = £ e >= a el a c ka + e dd + *£(V ka cl a d + V k * a ^ Cka ). ( 8 0) 
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In the equation-of-motion technique we need the following commutators: 

[Ck a , H] = ZkaCka + V^d 



The equation-of-motion for Pi is 



= -E^^+Ew^)- 

The hybrid functions have the following equations of motion: 

^< C L<*> = -«^{4^-^a{^) + E V *'a'< C L*a'> + eo{4 a d) 



or, in Fourier space, 



dt 

ifo{d*c ka ) = € ka {Sc ka ) + V ka {dU)-Y^V klal {cl a ,c ka )-s a {Sc ka ), 



[u-(e Q -e ka )}{cl a d) a - -V k ' a {dtd) u + V k * a 2n5(u>)f a ( eka ) 
[w-Cete-eoJJ^Cto}.* = Verity- - Vfc a 2^(w)/ a (6 te ). 



(81) 



(82) 
(83) 

(84) 
(85) 

(86) 
(87) 



Here we used the fact that the contacts are noninteracting and in equilibrium, in order to write {c\ a c k - a i) = 
Skk'S a ,a'fa(eka)- Substituting these in the Fourier-transform of (83) gives " 



w{Sd) a = -El^l 2 [Z*5(w)f a (e ka ) - (Sd) a ] 



ui + irj - (e - e ka ) w + it] - (e ka - e ) 



(88) 



where we by hand set u) w + it]. We try to proceed within the spirit of WBL. Consider first the term proportional 
to f a (eka) (call this (1)). Because of 6(u), the principal value contributions cancel and we find 

(1) -* -2*d(w) £ \V ka \ 2 f a (e ka )(-2m)6(e - e ka ) = 2m5(w) ^ a U(e ), (89) 

ka a 

where we identified T* = J2h ^\ v ka?5{€o - e ka ). The second term is evaluated in a similar fashion: 



\ui + irj - (e - e ka ) u + iT)-(e ka -e )] N ^ ' 



(90) 



where we again used the WBL to suppress the principal value parts, and to ignore the energy-dependence of the T's. 
Putting the pieces together, we have arrived at 



u>(d!d) u = iJ2T a P*S(«,)f a (e ) - (d*d) u }, 



(91) 



Going back to time-domain with / du/{2ir) exp(-uot), introducing the P's, and recalling that P 4- Pi = 1, allows us 
to write 



6Fl 



(92) 



in perfect agreement with the ME. 
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Derivation of the Landauer formula using Nonequilibrium Green functinos 

A. P. Jauho 



I. BACKGROUND 

During the lecture we discussed the Meir-Wingreen formula for calculating current through a mesoscopic region 
(where interactions may be present). The full details are given in Section VII in the handout note. Here we use a 
similar approach to rederive the famous Landauer formula. This is probably the most complicated and clumsy way 
of deriving a well-known result, nevertheless it is an illustrative example. The formula emerges "almost naturally" 
if we just remember to stick to the spirit of the scattering theory for conductance: the leads are described in terms 
of scattering states, i.e., free fermions. Your task is to go through the derivation sketched below, filling details, and 
discussing the physics. 

H. CALCULATION OF THE CURRENT FOR A TUNNELING STRUCTURE 

Just as in the general case, the current from the left contact through left barrier to the central region can be 
calculated from the time evolution of the occupation number operator of the left contact. Let us first introduce some 
notation. The contact Hamiltonian is H c = J? Left -f flujght, where 

L 

•Hkight = X) €flC fl CB » (1) 

R 

where we introduced quantum numbers L, R to label the states in the contacts (these quantum numbers contain all 
the necessary labels, such as longitudinal momentum, sub-band, or spin indices). The tunnel Hamiltonian H T is 

* H T = ^2[tLRc{cR + t RL c R c L ]. (2) 

RL J 

Task 1. Calculate the current operator 1-^ by forming the commutator [N Leit , #]. 
Task 2. Show that the expectation value for the current is 

JL S (/Left) = -i [tLR(c{c R ) - t RL {c R d,)\ • (3) 
RL J 

Next we need to define some Green functions (we set ft = 1, and reintroduce it in the final expression for the 
current): 

G RL (t-t') s i{cl(t')e R (t)) , 

Gi R {t-t) = i(c R (t')c L {t)) . (4) 

We see that the current is given by the time-diagonal components of the Green functions defined in (4). These 
functions have the property Gf R (t, t) = - [G^ L (t, t)]*. 



Task 3. Show that the current can be expressed as 

2e 

He] 

LR 



Jl = j Re{^ii^(«, i)} . (5) 



Next, one needs an expression for G^t - 1'). For the present case, with non-interacting leads, a general relation 
for the contour-ordered Green function Grl(t, t 7 ) can be derived rather easily with the equation-of-motion technique. 
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TABLE I: Rules for analytic continuation 



Contour Real axis 



c = 




C< = f t [A r B< + A < B a ] 






CT = J t A r B r 


D = 


Sc ABC 


D< = UA T B r C< + A r B<C a + A < B a C] 






D T = J t A T B r C T 


C(T, 


t') = A(t,t')B(t,t') 


C<(t,H)=A<(t,t')B<(t,1?) 






C r (t,t') = A<(t,f)B r (t,if) + A T (t,t')B<t,t,t') 






+A r (t,t')B T (t,t') 


D(r 


r')-^(r,T')B(r',r) 


D<(t,t?) = A<(t,lf)B>(1?,t) 






D<-(t,t') = A<(t,i?)B*(lf,t) + A r {t,t')B<{t',t) 



Since the nonequffibrium theory is structurally equivalent to equilibrium theory, it is sufficient to consider the T = 
equation-of-motion for the time-ordered Green function G^t ~if) = -i(T{c R t)c\{t' )}): 

Task 4. Show that the equation of motion is 

-*) = £i Gk(t - f) + £ GWt - H)t R , L , (6) 

R' 

where you find the right lead time-ordered Green function function G^i - *') = -i{T{c R (t)c a ,(t')}). 

Note that it is crucial that the leads be non-interacting; had we allowed interactions in the leads, the equation-of- 
motion technique would have generated higher order Green functions in (6), and we would not have a closed set of 
equations. In the spirit of the scattering theory of conductance, one now must use a noninteracting version of Grrk 

Grri — * SRRigR. 

Question. Without peeking further down the text, do you remember how the retarded noninteracting Green 
function g R looks like? 

The equation of motion now simplifies: 

JjGkt* - = e L GR L (t - f) + g R (t - f)t RL , ( 7 ) 

We can interpret the factors multiplying G t RL {t -If) as the inverse of the (left) contact Green function operator, and 
introduce a short-hand notation: G^g^ = g%t RL . By operating with g\ from right, we arrive at 

Gjidt - = / dt l9R (t - tijtjutflfo - f ) . (8) 

Since the equilibrium and nonequilibrium theories are topological^ equivalent, this equation has in nonequilibrium 
precisely the same form, except that the intermediate time integration runs on the complex contour: 

G rl (t,t') = j dT 1 g R (T,T 1 )t* RL g L (T 1 ,T') .' (9) 

The analytic continuation rules we need are listed in Table 1 and can be readily applied to complete the next task. 

Task 5. Show that 

G RL (t -t') = J dht RL \g R {t - tdgHti. - 1') +g&t- t{)gl{h - 1')} , (10) 
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where the Green functions g^'^ for the leads are familiar: 

ai/ntt ~ = i{<i /R ^)c L / R (t)} = ih/ R (£L/ R ) exp [ - is L/R (t - 1')] , (11) 
- = =T»(±* =F *0«Ci/«(*). «U(*0}) = =FW(±t =F *') exp [ - fe^t - i')] • (12) 

Here }l/r{£l/r) = [exp[(ex,/jj - Ht,/ R )/k B T] + 1] -1 is the equilibrium distribution in a given lead. 
The Fourier transform of (10) is 

G RL {e) = t RL \g R {e)gi{e) + gfc)gi(e)] . (13) 

Task 6. Go through the details of the following calculation. 

J i = |E/ %M*L&L\f&)ti{e) + gZ(e)fi(e)]} 

= |E *«*Iui / * [fr - /*(«*)] *( e - e ^( e - «*) 

= YE^^t^( e i)-/«( £ fl)] 5 ( e H-et)- (14) 
Finally, we evaluate the sums over L/R by making the replacements 

Yl-'ffeL/R E ' t LIR-*'tn L ln R { s L/R), (15) 

where n£/# contains the discrete quantum numbers. Using these, you can complete the last task: 

Task 7. Show that the current becomes 

J = j J (feTr[ttt(e)][ni(e)-na(e)l 

= x Tr[ttt( ej? )]y, (16) 

where the second equality holds in linear response, and at low temperatures. This is, of course, the celebrated 
conductance formula, 

III. ADDITIONAL PROBLEMS 

Time permitting, here are additional small problems you may want to discuss during group work. 

• Examine the expression (71) in the hand-out note. In particular, take T R (e) = T L (e) = r (independent of 
energy), and A=0. Evaluate the current at T = 0, and discuss the rather boring result. What should one do to 
obtain a true resonant behavior, as seen in the experiments? 

• Examine the Master Equation calculation of the resonant tunnelling current (ending in Eq.(79)). . 

• Go critically through the NEGF derivation of the ME (p.13 in the note). 



